NEW USES FOR SENSITIVITY ANALYSIS: HOW DIFFERENT MOVEMENT 

TASKS EFFECT LIMB MODEL PARAMETER SENSITIVITY 


INTRODUCTION 


Sensitivity analysis techniques have traditionally been used by the 
systems engineer to help understand the behavior of complex systems. 
Although the detailsof the techniques seen in practise differ, the basic 
approach is the same: A system parameter is varied in a controlled, 
systematic manner, and the subsequent variations in output are measured and 
described quantitatively (Frank (1978), Tomovic and Vukobratovic (1972), 
Lehman and Stark (1982)). When the system is nonlinear, as is the usual 
case, numerical techniques involving computer simulation typically need to 
be employed. Insights gained from such techniques are of value both for 
systems analysis and design. 

The methods presented here represent extensions of previous work (Clark 
and Stark (1976), Hsu et al (1976), Lehman and Stark (1979), Bahill (1980) 
and Zangemeister et al (1981)), only with a wider range of sensitivity tools 
employed. Furthermore, the model considered here is for limb flexion- 
extension movements, rather than for head or eye rotation. The model 
structure has also been expanded and the constitutive equations representing 
basic muscle properties improved so as to more accurately characterize basic 
neuromuscular system dynamics. Consequently, there are a larger number of 
internal model parameters. A greater number of output behaviors are also 
considered. 

In addition to presenting this expansion of previous sensitivity 
analysis tools and extending these methods to a larger number of model 
parameters, a major role of this presentation is to show how sensitivity 
analysis results are a function of the model task. When the task under 
consideration is changed (i.e. the model input controller signal sequence is 
fundamentally different, resulting in a different type of output), the 
relative role of each parameter in affecting performance also changes. This 
fact, surprisingly neglected in the literature, is developed quantitatively 
here. 


The result is one model that can adequately simulate any basic 
physiologically realizable flexion-extension task and a set of sensitivity 
tools that help explain the relative role of any specific parameter for any 
particular task - tools that can help make the goal these modeling efforts, 
gaining insight into the role of biomechanical systems in neuromotor 
control, a reality. 
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METHODS: 


A. Model Structure: Before presenting the sensitivity analysis 
protocol employed here, it. will first be advantagous to develope a basic 
understanding of the biomechanical system being modelled. The first step to 
any modeling effort is to assume a basic structure. Once chosen, this 
structure will be the major constraint on the success of the model - too 
simple a structure can result in a poor approximation of actual behavior, 
while too complex a structure reduces Insight (or results in a model with 
poorly defined parameters). The basic structure for the model, based on the 
classical muscle work of Hill (1938) and elbow flexion-extension work of 
Wilkie (1950) and supported by numerous more recent experimental work on 
muscle mechanics, is presented in Figure 1. The sixth-order strucuture has 
been found to be the lowest order structure that is capable of approximating 
all fundamental muscle properties needed for an antagonistic pair of lumped 
"equivalent" muscle actuators rotating a joint. 

FLEXOR NEURAL INPUT EXTENSOR NEURAL INPUT 


/ 

/ 

/ 

/ 

/ 

/ 

/ 

/ 

/ 

/ 

/ 

/ 



« 

® 


® 

® 


FIGURE 1: Model of System Showing the Nonlinear Blocks. 

Lumped Flexor Muscle is on the left. Extensor on right. 

Experimental work is often able to approximately isolate each of these 
elements in the model. Fundamental to such an approach is this idea of an 
"equivalent muscle". This concept of an "equivalent" muscle for the lumping 
of a number of synergistic muscles has been previously developed, based on 
experimental work, for both the flexor group (Bouisset et al (1973, 1976) 
and the extensor group (Cnockaert and Pertouzon (1974)) and confirmed by 
Cnockaert (1978) and Le Bozec (1980). This idea is supported and expanded 
on here in the following sense: not only is it a good representation for 

the ideal case of elbow flexion-extension but it also should be expected to 
hold for more complex one degree of freedom movements such as wrist or head 
rotation because two lumped antagonistic muscles with the blocks described 
above should be structurally capable of approximating all basic muscle 
properties for such movements whenever muscles contract approximately 
synerglstically. In these more Involved movement systems, however, parameter 
identification is more difficult. 

A good summary of much of the work on the material properties of muscle 
is found in the review by Close (1972). Once fiber type and fiber 
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orientation are determined, basic skeletal muscle material properties for a 
given muscle can be well estimated. Geometric data for the muscles around 
the elbow joint also exists (An et al (1981), Amis et al (1979)). By 
creating algorithms that combine material and geometrical information, first 
generation parameter values for the torque-velocity and series elastic 
elements can be established. These results are then combined with the 
wealth of experimental work on controlled Intact limb movements, the best of 
which include Bern et al (1947), Wilkie (1950), Pertuzon and Boulsset 
(1973), Jorgensen et al (1971, 1976), Hatze (1981a, b) and Komi (1973) for 
torque-velocity information, Wilkie (1950), Goubel and Pertuzon (1973) and 
Cnockaert et al (1978) for the series elastic relation, and Boon (1973) and 
Hayes and Hatze (1977) for passive viscoelastic data. Limb inertial data 
exists in abundance. Insights from preliminary sensitivity analysis work 
(not presented here) are also used for fine-tuning parameters. The actual 
protocol followed for parameter development is beyond the scope of the 
present presentation and will not be described here. The resulting model 
parameter values for the elbow flexion / extension model, one of the five 
models currently under pursuit, are displayed in Table 1. 


PARAMETER: VALUE: CONSTITUTIVE EQUATION: 

Passive Plant: 

Jp: 0.06 Kg-m**2/rad (or N-m-sec**2/rad) 

Bp: 0 .15 N-m-sec/ rad 

Kp : 1.4 N-m/ rad 

Kpl: 0.0001 ... 

Kp2: 10.0 

Series Elasticity: 

Ksl-f; 4.8 N-m/ rad 

Ks2-f: 7.0 

Ksl-e: 4.5 N-m/rad f Fks = Ks. *(exp(Ks 2 *(xh-x))-l) 

Ks2-e: 7.2 ' 

Torque-Veocity : 

Af-f: 0.34 

Bh-f : 8.0 rad/ sec 

Af-e: 0.30 

Bh-e : 7.0 rad/ sec 

Pm-fv: 0.3 

Af-fv: 0.6 

Bh-fv: 0.25 

Activation Dynamics: 

Tal : 40 ms 

Ta2: 10 ms 

Fmax-f : 60 N-m 

Fmax-e: 50 N-m 



^ Fkp = Kp*x + Kpj^*(exp(Kp 2 *x)-l) 


TABT.K 1: Current Parameter Values for the Elbow Flexion-Extension Model. 
Constitutive Equations are for: Parallel Elasticity (Fkp), 

Series Elasticity (Fks (=Fm)) and Torque-Velocity (Bm). 
(Parameter values are for 70 Kg male of average strength.) 
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The passive plant values represent the inertia and viscoelastic 
properties of the lumped joint/muscle system, with all three elements in 
parallel as is usual. The exponential fit for the parallel and series 
elasticities has become a standard representation for load-bearing 
collagen-elastin based soft tissue (Fung (1968), Glantz (1974) and Hatze 
(1981)). The "Hlll°s” parameters (Af and Bh) for shortening muscle are the 
standard representation for the classic force-velocity relation of muscle, 
used constantly in the literature to document experimental results (see, for 
instance. Close (1972)). The scaling of the Instantaneous torque-velocity 
relation by the activation level was suggested as early as 1956 by Wilkie, 
and has been supported by the work of Pertuzon and Bouisset (1973). The 
relation used for lengthening muscle (an Inverted, skewed Hill°s-^type 
equation) is a new method that appears to adequately approximate past data 
(Joyce et al (1969a, b), Komi (1971) and Hatze (1981)), plus numerous 
observations that the peak eccentric torque is about 30% above the peak 
Isometric force. Activation dynamics is simulated by two time constants, 
compatible with the basic neuromuscular literature (see Close (1972) or 
Bahler (1967) for reviews). This second-order form represents 
simplifications suggested by the more detailed work of Lehman (1983) and 
Hatze (1981). There is also numerous isometric peak torque data available 
in the physical education literature - the values presented here are for a 
"typical" human male. For reasons of clarity, the static torque-angle 
parameters (based on an abundance of literature) were not presented above. 

As seen above, all indications are that all of these elemental building 
blocks are nonlinear. The function of these nonlinear properties is still 
poorly understood, and one of the main problems faced is to explore the 
sensitivity of the system to these nonlinearities. There is ample evidence, 
supported here, that the relative importance of various parameters is a 
function of the task in question. Consequently, it can be a major mistake 
to over-simplify this basic system if one is interested in a variety of 
movement tasks. Furthermore, since sensitivity methods provide just the 
information needed for task-specific model simplification, it is suggested 
that the more complex model be considered first - any model simplification 
is then based on a solid foundation. 

B. Computer Simulation Algorithm: The simulation algorithm is 
contained within a more general set of modules that are linked to a main 
routine, called "JAMM" (Juiced-up Antagonistic Muscle Model). This user- 
friendly program will simulate second, sixth and eighth order models with 
degrees of nonlinearity ranging from linear to highly nonlinear. Once the 
biomechanical model of interest is chosen, a data base, complete with all 
the current numerical values of parameters for any user-desired combination 
of linear/nonlinear parameter defaults, is read. The user is prompted for 
parameter modification, for various external loading options, for the type 
of run (interactive, sensitivity analysis, optimization), and for the 
controller signal input sequence for each equivalent muscle. 

The options under sensitivity analysis include; determining the 
parameters that are to be varied, one by one, for a given run and 
determining the range of the parameter variation and the number of times 
varied. Parameter variation is by a reciprocal format (for example, 4/5 and 
5/4 of nominal). Typically results for five reciprocal pairs are obtained. 
Raw behavior data, behavior and parameter ratios, and linear and logarithmic 
sensitivity coefficients are stored for later plotting and/or printing. 
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Input Sequence 


BI(»fECHANICAL SYSTEM 


Joint Output 


— >• Behaviors 


4 Parameters to be Varied ^ • S^^^ =» (dB/Bo)^/(dP/Po) ^ 


Figure 2 Schematic of Sensitivity Analysis Method. 

Note that varied parameters can include system, input or 
disturbance values, and that input sequence depends on task. 


C. Sensitlvty Analysis Protocol: Past sensitivity analysis work on 
eye and head systems concentrated on the development of a "sensitivity 
matrix" (Hsu et al (1976), Lehman and Stark (1979), and Zangemelster et al 
(1981)). A schematic of this basic method is presented in Figure 2. An 
input sequence that will define a certain type of task is chozen, and 
output trajecties are measured. In general, the input and output can be 
scalar or vector quantities. Here, for the generalized equivalent flexor 
and extensor muscles, there are two inputs, one to each of the lumped 
equivalent muscles. The outputs of Interest are the position, velocity and 
acceleration of the limb, plus the muscle torques. The "behaviors" of 
interest are a function of these output trajectories. The parameter varied 
is typically an Internal model parameter value, but may also be an input 
signal parameter value, such as a pulse height or width, or an external 
disturbance. Traditionally, each column would indicate the sensitivity of 
all the different behaviors (each on a different row) to that column'’s 
parameter. 


The actual value of each matrix element is called the "sensitivity 
coefficient", of the i-th behavior to the j-th parameter. This 
coefficient represents the relative change in behavior divided by the 
relative change in parameter (dB/B^)^/(dP/P ),, where Po and Bo are the 
nominal parameter value and the resulting nominal behavior value, 
respectively. The range of the change in parameter for the determination of 
the matrix coefficients is up to the discretion of the user. Typically, the 
range chosen for the sensitivity matrix coefficient computation was from one 


half to twice the nominal parameter value. Another design consideration is 
the equation used to determine the coefflcent. Two equations are used here: 


"linear": ((B2"B^ ) 

"logarlthittlc": (log(B2/B„)/log(P2/Po))- (log(Bj/B^)/log(P^/PQ)) 


where P 2 is the parameter value greater than nominal and B 2 the resulting 
behavior value; P^ the reciprocal fraction of P 2 and Bj^ the resulting 
behavior due to Pj^. Because the first method gives a value proportional to 
the relative difference in behavior without regard to one direction maybe 
having a greater shift, it tends to weigh behavior changes greater than B 
disproportionally more than those below. The second method weighs ratios 
both below and above nominal equally. For this reason, the second method 
is usually prefered. Notice that, if the behavior where to change in a 
manner proportional to the parameter change, the sensitivity coefficient for 
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either method would be "1.0". 

Once this "sensitivity matrix" is completed, it gives a global view of 
model behavior for the task under question. Table 2 presents the results for 
a simple "generic" run. Here, an input signal, about 30% of maximum, is 
applied to the flexor group for 200 ms. The extensor group is about 3% of 
maximum. The model is run for each of the parameters chosen for variation at 
values one half and twice nominal. The behaviors of Interest are measured, 
and the sensitivity coefficients determined, here by both methods. 


Each of the resulting column gives one a feel for how a given parameter 
effects the various behaviors, while a given row indicates what parameter(s) 
most influence the particular behavior. For convenience in matrix 
inspection, the following conventions are used; the highest value in each 
column is printed in italics; the highest in each row is in boldface; and 
the three three most influencial parameters are also printed in boldfaced 
italics. For this example coefficients for both the "linear" and 
"logarithmic" decriptions are provided. Note the similarity in coefficient 
values. All later work uses only the logarithmic method of determination. 


TABLE 2: 


SENSITIVITY MATRIX FOR TASK; “Generic, Medlia-Speed Hovemenc"; 



NOMINAL 

Jp Bp Kp 

Kpl 


Ksl-{ 

Ks2-f 

Af-f 

Ib-f 

Fvmax 


Magnj 

114 deg 

-0.012 -0.108 -0.026 
-0.013 -0.120 -0.029 

-0.012 

-0.013 

- 0.289 

- 0.398 

-0,002 

-0.002 

-0.004 

-0.004 

-0.200 

-0,223 

0.037 

0.055 

-0.101 

-0.112 

-0.046 

-0.050 

Vnax: 

565 d/s 

-0.141 -0.125 0.005 

-0.162 -0.138 0.006 

-0.000 

-0.000 

- 0.367 

- 0.575 

-0.007 

-0.070 

-0.019 

-0,020 

-0.290 

-0.333 

0.542 

0.593 

-0.059 

-0.064 

-0.115 

-0.131 

Amax: 

5492 d/s/s 

-0.510 -0.057 0.077 

-0.556 -0.062 0.082 

-0.002 

-0.002 

1.732 

0.924 

-0.021 

-0.023 

-0,039 

-0.043 

-0.096 

-0,105 

0.274 

0.300 

-0.057 

-0.063 

- 0.415 

- 0.468 

Amin: 

-4858 d/s/s 

- 0.671 - 0.231 0.150 

- 0.757 - 0.227 0.156 

-0.028 

-0.030 

1.960 

1.000 

-0.020 

-0.021 

-0,041 

-0.044 

- 1.128 

- 0.800 

5.573 

1.749 

-0.002 

-0.002 

-0.390 

-0.454 

Fm-f: 

9.1 Nm 

0.154 0.005 -0.015 

0.166 0.006 -0.017 

-0.001 

-0.001 

0.295 

0.264 

0.005 

0.006 

0.022 

0.040 

-0.063 

-0.069 

0.185 

0.202 

0.033 

0.050 

-0.123 

-0.134 

Fb-f: 

14.4 N-m 

-0.100 -0.056 0.011 

-0.110 -0.061 0.012 

-0.000 

-0.000 

- 0.364 

- 0.569 

-0.001 

-0.020 

-0.004 

-0.005 

0.106 

0.114 

-0,170 

-0.191 

-0.029 

-0.032 

-0.138 

-0.160 

Tmagn: 

363 ms 

0.298 0.097 -0.064 

0.302 0.104 -0.071 

- 0.041 

- 0.044 

0.000 

0.000 

-0.009 

-0.010 

0.017 

0.080 

0.000 

0.001 

0.233 

0.241 

0.017 

0.018 

0.252 

0.252 

Tvmax: 

209 ms 

0.131 -0.013 -0.035 
0.147 -0.014 -0.038 

0.000 

0.001 

0.077 

0.078 

0.000 

0.001 

0.003 

0.003 

-0.061 

-0.068 

0.042 

0.045 

-0.003 

-0.003 

0.057 

0.061 

Tamax: 

59 ms 

0.260 -0.023 -0.011 
0.278 -0.024 -0.012 

-0.011 

-0.012 

-0.158 

-0.192 

-0.068 

-0.073 

-0.192 

0.207 

-0.034 

-0.036 

0.068 

0.071 

0.045 

0.048 

0.294 

0.320 

Tamln: 

278 ms 

0.103 -0.132 -0.007 
0.111 -0.135 -0.080 

-0.031 

-0.033 

- 0.528 

- 1.182 

-0.005 

-0.005 

-0.017 

-0,018 

-0.120 

-0.125 

0.029 

0.032 

- 0.134 

- 0.136 

0.089 

0.096 

T£m-f: 

54 ms 

0.284 0.025 -0.037 

0.390 0.027 -0.040 

0.000 

0.010 

0.877 

0.605 

- 0.086 

- 0.095 

- 0.247 

- 0.260 

-0,037 

-0.040 

0,099 

0.105 

0.037 

0.040 

0.247 

0.265 

Tfb-f: 

273 ms 

0.023 0.000 -0.011 

0.050 0.000 -0.012 

0.000 

0.010 

0.010 

0.011 

0.000 

0.010 

0.000 

0.001 

0.000 

0.000 

0.000 

0.000 

0.000 

0.001 

0.036 

0.040 
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Use of a coefficient determined by two values only gives the linear 
slope over the operating range of the two parameters. For a nonlinear 
system, this may give misleading information (discussed later). For this 
reason, the parameter ratio range for a given sensitivity matrix is an 
Important design variable. Thus, if one is interested in a deeper 
understanding of the role of a certain parameter, a simple column of 
coefficients is not enough. 

The tools desribed here for a more in-depth examination of the role of 
a specific parameter will be called "sensitivity graphs" and "sensitivity 
trajectories". They are best used as the next step after a preliminary 
sensitivity matrix has been developed for the task in question. The 
"sensitivity trajectory" is simply the set of output versus time plots that 
result from a range of parameter variation, in superimposed plots (Figure 3, 
left panel). Inspection of these output plots can be a surprisingly 
effective way of coming to an understanding of the role of the parameter, 
making use of human talents for visualizing information, and putting the 
column of sensitivity coefficients in proper perspective (these 
coefficients can be occasionally misleading (discussed later). This simple 
step should be used on all parameters with significant sensitivity columns. 

"Sensitivity graphs" further expand ones Insight into model sensitivity 
to a certain parameter of interest, and also bring together sensitivity 
columns and trajectory information. This method consists of graphing 
behavior ratios versus parameter ratios for a wide range, as in Figure 3 
(right panel). Notice that each "graph" is basically a graphical extension 
of each coefficient, basically showing the five possible coefficients 
(slopes) that could be placed in the particular location. Typically 
logarithmic scales are employed. Visual inspection of this graph provides 
information on how linearly the behavior changes with parameter variation. 
It also suggests the useful operating range of the parameter of Interest. 
This is possibly the most important sensitivity tool from a design 
perspective. 
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Figure 3: Sensitivity Trajectories (left panel) and Sensitivity Graphs 

(right panel) for the controller parameter PH-f, the piilse 
height for the agonist pulse. Task #1. Large dash is for 4/5 
and 5/4 times nominal. For Sensitivity Graphs, range is l/lO 
to 10, on logarithmic scales, for both the parameter (abscissa) 
and the behavior (ordinate). 
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RESULTS: 


Task #1: "Generic" Run: Some simple results, for a "typical", moderate- 
speed run, were displayed In Table 2 (sensitivity matrix) and Figure 3 
(sensitivity trajectories and sensitivity graphs). From the sensitivity 
matrix, it is seen that the most important system parameters are the plant 
inertia (Jp), especially for peak accelerations and for the timing of all 
peak values; Hill°s parameter Bh-f, especially for peak magnitudes; the 
time constant of activation (Taj^), and the parallel elastic fit parameter 
Kp 2 - From the nature of the equation for the parallel elasticity (Table 1) 
we see that the location of the parallel elastic concavity will 
automatically define the position operating range. The model output is 
also particularily sensitive to controller signal pulse parameters such as 
PHI and PH2 - not an unexpected finding since a well-designed tracking 
system is usually sensitive to its own input for tracking tasks. 

Based on these sensitivity matrix results, the parameters mentioned 
above appear to be of particular Interest for this bask. In Figure 3, 
variation in the agonist pulse height parameter (PHI) is displayed using 
sensitivity trajectories and sensitivity graphs. These results show more 
explicitly the effect of varying the agonist pulse height. For reference, 
an average adult male can contract the flexor group to about 60 N-m. 
Similar plots, not presented here, are then produced for the other highly 
sensitive parameters. 
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Task #2; Unloaded Fast Voluntary Elbow Flexion: The effects of a 60 N- 
m, 100 ms agonist pulse are displayed in Figure 4. The corresponding 
antagonist pulse during this time is only 1 N-m., Notice that the peak 
angular velocity is about 14 rad/sec. (For a longer pulse width, one finds 
a peak velocity of about 19 rad/sec both for the model and for the average 
intact adult human (Bern (1947), Wilkie (1950), Pertuzon and Bouisset 
(1973)). Zero degrees is defined here as the rest elbow position (where 
the lumped parallel elasticity is zero), of 100 degrees between the humerous 
and the ulna. Notice that movements of about 80 deg are possible in each 
direction, with the nonlinear parallel elastic element automatically keeping 
the joint position within this physiological operating range. 


In Table 3 a more complete sensitivity matrix for Internal system 
parameters. Including flexor and extensor parameters for both series 
elastic and torque-velocity (shortening and lengthening) properties, is 
presented. Notice that, while the general trends in this table are similar 
to those in the previous table, the details are quite different. Also 
notice the general insensitivity of the "extensor" muscle parameters. This 
shows that it is the flexor parameters that are of primary importance for 
this particular task and furthermore shows the relative significance of the 
various flexor parameters to each other for each behavior. 

In Figure 5 sensitivity trajectories and sensitivity graphs are 
obtained for a few of the more important internal system parameters such as 
Bh, Jp and Tal. Notice that, by combining all three techniques, a 
remarkably clear picture of the role of each of these parameters emerges. 
We see, by all three methods, that the inertial term affects mainly 
acceleration information and the timing of peak values. The, latter two 
methods both show that the system is more sensitive to increases in inertia 
- something common in everyday life and sporting events. All three methods 
also show that the torque-velocity parameter Bh (Figure 5b) mainly 
influences magnitude information, with less effect on timing. The 
acrtivation time constant parameter Ta^, which basically filters the 
neuromuscular signal before "passing" it, effects system behavior as might 
be expected. Notice that the sensitivity increases proportionally more when 
the parameter increases in value than when it decreases (best seen by the 
sensitivity graph). 


TABLE 3: 


SENSlTIVm lunu for task: Unloaded Movaatnt*: 



ROHIKAL 

Jr 

Bp 

KP 

Kpl 

Kpl 

Kal-t 

Ka2-f 

Kal-a 

Ka2-a 

Af"f 

Af-a 


Bb-a 

ta-fa 

Af-fv 

Bh-fv 

141 

Ta2 

Mago: 

122.6 dag 

0.050 

-0.065 

-0.036 

-0.015 

-0.200 

-0.002 

-0.006 

0.001 

0.002 

-0.173 

-0.017 

0.335 

0.023 

-0.073 

-0.017 

0.022 

-0.001 

-0.003 

Vuz; 

805 d/a 

-0.247 

-0.090 

0.024 

-0.004 

-0.130 

-0.017 

-0.055 

-0.000 

-0.003 

-0.323 

-0.005 

0.565 

0.005 

-0.031 

-0.005 

0.0U5 

-0.226 

-0.038 

Aaax: 

16A07 d/a/a 

-0.527 

-0.051 

' 0.037 

-0.004 

-0,170 

-0.015 

-0.011 

-0.005 

-O.OU 

-0.166 

-0.004 

0.387 

-0.008 

-O.D26 

-0.004 

-0.008 

-0.447 

-0.112 

Mint 

6081 d/a/a 

-0.766 

-0.005 

0.080 

-0.002 

-0.046 

-0.037 

-0.095 

0.000 

0.000 

-0.905 

0.022 

1.932 

-0.026 

0.137 

0.022 

-0.026 

-0.528 

-0.051 

Taagn: 

320 •• 

0.131 

0.750 

-0.101 

•^.069 

-0.755 

0.022 

0.057 

0.000 

0.000 

0.000 

0.000 

0.087 

0.000 

-0.017 

0.000 

0.000 

0.333 

0.049 

Tvux: 

125 u 

0.160 

-0.029 

-0.015 

0.000 

0.036 

-0.012 

-0.029 

0.002 

0.009 

-0.094 

0.022 

-0.136 

0.002 

0.770 

0.022 

0.023 

0.137 

•0.038 

Taaax: 

51 «a 

0.281 

-0.028 

0.000 

0.000 

0.170 

-0.015 

-0.287 

0.000 

0.000 

0.000 

0.000 

-0.014 

0.000 

0.042 

0.000 

0.000 

0.252 

-0.112 

Taain: 

189 M 

0.234 

-0.061 

0.006 

0.0038 

0.078 

-0.004 

-0.011 

0.000 

0.000 

H).284 

-0.008 

0.134 

0.000 

0.038 

-0.008 

0.000 

0.190 

-O.OU 

tw-fi 

15.6 

0.283 

0.009 

-0.013 

-0.004 

-0.069 

0.003 

0.040 

-0.004 

-0.010 

-0.145 

0.006 

0.343 

-0.017 

0.036 

0.006 

-0.017 

-0.213 

-0.060 

Fv*a: 

8.4 N-a 

0.080 

0.009 

0.078 

0.019 

1.56f 

0.000 

-0.002 

-0.000 

-0.010 

0.026 

0.070 

-0.032 

-0.043 

0.183 

0.037 

-0.043 

-0.070 

0.002 

Pb-f: 

37.6 N-s 

-0.219 

-0.028 

0.000 

-0.001 

-0.063 

0.016 

0.059 

0.001 

0.001 

0.101 

-0.003 

-0.136 

0.004 

-0.019 

-0.030 

0.004 

-0.424 

-0.108 

Fb-a: 

1.9 N-b 

0.101 

-0.027 

-0.008 

0.001 

0.741 

-0.004 

-O.OU 

0^000 

0.010 

-0.338 

0.156 

0.712 

-0.167 

0.971 

0.156 

-0.167 

-0.118 

-0.010 

Fa-fa: 

11.4 h-m 

0.500 

0.018 

-0.086 

-0.034 

-0.110 

-0.029 

-0.057 

-0.005 

-O.Olt 

-0.202 

-0.040 

0.468 

-0.007 

-0.027 

-0.004 

-0.008 

-0.488 

-0.U2 

Tf«-ft 

50 M 

0.287 

0.015 

-0.015 

0.000 

0.000 

-•.102 

-0.310 

0.000 

0.000 

0.000 

0.000 

0.014 

0.000 

0.029 

0.000 

0.000 

0.240 

0.171 

Tfa-e> 

207 Ba 

- 

0.028 

- 

- 

- 

0.003 

0.003 

6.000 

-0.003 

0.084 

0.000 

• 

-0.021 

- 

0.000 

-0.021 

0*279 

0.082 

Tfb-fs 

104 M 

0.007 

-0.007 

0.000 

0.000 

-0.022 

•0.007 

-0.007 

0.000 

0.000 

0.007 

0.000 

0.000 

0.000 

0.000 

0.000 

0.000 

0.041 

0.041 

Tfb-a: 

193 «a 

O.lSl 

-0.008 

-0.019 

0.000 

-4.492 

0.000 

-0.004 

0.000 

-0.004 

-0.328 

-0.004 

0.159 

-0.049 

-0.494 

-0.004 

-0.049 

0.215 

0.077 

W-fa: 

55 Ba 

0.249 

0.026 

0.000 

0.013 

“ 

-0.080 

-0.292 

0.001 

0.000 

-0.026 

0.000 

0.013 

-0.013 

0.039 

0.000 

-0.013 

0.287 

0.193 
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Figure h: Plots of input, output, and internal variables for the 

"normal" r\m of Task §2. Lower left contains the model 
input variables (FK-f & FN-e) and the resulting muscle 
output torques (FM-f & FM-e). Upper and middle left contain 
kinematic output infonaation as well as the positions and 
velocities of flexor and extensor internal nodes (dashed). 

On the' right, from top to bottom, are internal variable 
plots for the series elastic element, the instantanous 
externally seen torque-velocity behavior, the viscous 
muscle torque versus node velocity, and the torque propogation 
for varous model elements. For the 3 top right plots, time 
is an implicit parameter, with a point being produces every 
10 ms . 
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Figure 5aJ Sensitivity Trajectories (left panel) and Sensitivity 
Graphs (right panel) for parameter Jp (plant inertia) , 
Task #2. Large dash is for 4/5 and 5/4 times nominal, 
smaller dash is for 1/5 and 5 times nominal. For 
Sensitivity Graphs, range is 1/10 to 10, on logarithmic 
scales, for both the parameter (abscissa) and behavior 
(ordinate) . 
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TORQUE <N-.) ACCEL (rod/«c/.«=) VELOCITY (rod/..e) POSITION (d.g> 


100 | 




Ttm-f X F0-m « 


Fn-f. Fm-f 




Fn-«, F«-b 


TIME (MC) 


Sensitivity Trajectories (left panel) and Sensitivity 
Graphs (right panel) for parameter Bh (one of the 
’’Hill’s" constants for shortening muscle). 

Large trajectory dashes: 4/5 and 5/4 of nominal, 
smaller dashes: 1/5 to 5 times nominal. 

Sensitivity Graph range: 1/10 to 10, logarithmic units 
for both parameter (abscissa) and behavior (ordinate) . 




TOmUE (H-iO ACCEL <r^amc/ame> VELOCITT Crod/Mc) 


TIME (Me) 



sensitivity Trajectories (left 

Graphs (right panel) for parmeter Tal " 

two time constants representing 

Large trajectory dashes. an consitivitv 

smaller dashes: 1/5 and 5 times 
Graph range: I/IO to 10, s' 

parameter (abscissa) and behavior (ordinate). 
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Task #3: Maximum Voluntary Flexion with an External Load: This task 
represents an extension of the previous task, only with the addition of an 
external "isotonic" load of 12 N-m and a 0.1 kg-m**2 added inertia (a light 
bar being grasped). The task is based on some of the classical experimental 
data of Dern et al (1947), and the results are consistent with this work. 
The time of contraction is now 250 msec. Run dynamics are displayed in 
Figure 6a. An inspection of Table 4 shows that, once again, while many of 
the parameters are very similar, others differ sharply. Such observations 
guide further analysis of this task (not presented here). 

In Figure 6b we see the sensitivity to a new type of parameter - a 
disturbance (i.e. the external load). Here the "base" external load of 12 
N-m is varied as any other parameter would be. The sensitivity of the 
system to such an external load becomes well understood via inspection when 
all three methods are used. Notice that the factor of ten parameter change 
(to 120 N-m) is not attempted since it would injure a normal limb. 

In Figure 6c and 6d we see the sensitivity to another new type of 
parameter: an initial condition. In this instance it is an initial agalnst- 
movement velocity of 2 and 2- rad/sec. Notice that, in the absence of any 
initial inertial dynamics, the lengthening torque-velocity properties are 
able to easily compensate for the initial condition of Figure 6d, with the 
effect being negligible by the time of peak velocity, position or negative 
acceleration - only the peak positive acceleration is significantly 
affected, and then only for initial velocities of greater than 4 rad/sec. 
The sensitivity for the "with" initial velocity is also compensated for 
fairly well. These findings are in contrast to previous results for the 
fast, low inertia eye movement system (Winters et al (1983). 


TABLE 4; 


SEKSIIIVtiy HATXU tOS: Maxlul Coacractlon Agalnac Extarnal Load 



NOMINAL 

Jp Bp 

Kp 

Epl 

Kal-£ Ka2-f 

Af-f 

8k<f 

Fafv 

Tel 

Te2 

PH-f 

PU-e 

nr LOAD 

Vo- 

Vo+ 

tUgn: 

114 deg 

O.Oll -0.032 

-0.013 

-0.398 

-0.001 -0.006 

0.095 

-0.268 

-0.031 

-0.018 

-0.003 

0.135 

-0.175 

-0.365 

-0.021 

0.021 

VB8X1 

565 d/e 

-0.131 -0.059 

-0.000 

-0.575 

-0.003 -0.017 

-0.433 

0.788 

-0.015 

-0.096 

-0.016 

0.536 

-0.259 

-0.401 

-0.007 

0.004 

Am«X! 

5492 d/a/e 

-0.574 -0.024 

0.002 

0.924 

-0.006 -0.009 

-0.129 

0.368 

-0.101 

-0.419 

-0.081 

0.759 

-0.184 

-0.475 

0.306 

-0.210 

Aaisi 

-4858 d/e/a 

0.228 -0.183 

0.028 

1.000 

-0.299 -1.201 

-0.512 

1.839 

-0.907 

1.480 

0.159 

1.814 

-0.539 

-1.507 

-0.068 

0.384 

Tugoi 

363 u 

0.268 0.056 

-0.040 

0.000 

-0.005 0.019 

0.000 

-0.065 

0.024 

0.184 

0.025 

0.088 

0.740 

-0.196 

0.038 

-0.057 

Tvmx: 

209 aa 

0.066 -0.003 

0.000 

0.078 

0.000 0.028 

-0.031 

0.017 

-0.003 

0.052 

0.039 

-0.021 

-0.010 

-0.006 

0.003 

-0.006 


59 aa 

0.293 -0.010 

-0.011 

-0.192 

-0.057 -0.274 

-0.011 

0.022 

0.023 

0.269 

0.154 

-0.123 

-0.012 

-0.069 

-0.180 

0.157 

Twloi 

278 aa 

0.285 -O.lll 

0.033 


-0.020 -0.007 

-0.061 

-0.051 

-0.113 

0.141 

0.031 

-0.014 

-0.108 

-0.037 

-0.097 

-0.040 

Fa-fi 

9.1 N-B 

0.157 0.002 

-0.001 

0.264 

0.003 0.023 

-0.064 

0.294 

0.017 

-0.123 

-0.028 

0.536 

0.036 

0.170 

0.133 

-0.086 

Fr«t 

7.1 N-a 

0.029 0.000 

0.019 

3.388 

O.OOl 0.000 

0.006 

0.060 

0.079 

0.000 

0.000 

-0.043 

0.058 

0.000 

0.000 

0.112 

rb-£: 

14.4 M-a 

-0.092 -0.026 

-0.000 

-0.569 

-0.000 -0.004 

0.066 

-0.120 

-0.007 

-0.120 

-0.020 

1.260 

-0.122 

-0,317 

0.009 

0.000 

Fb'as 

1.5 N-a 

0.231 -0.048 

-0.001 

3.269 

-0.002 -0.018 

-0.051 

0.838 

0.392 

-0.098 

0.070 

0.542 

-0.003 

-0.406 

0.020 

0.314 

Fw-fat 

5.0 N-a 

0.229 0.004 

-0.012 

3.614 

-0.005 -0.018 

-0.075 

0.377 

-0.006 

-0.224 

-0.041 

1.220 

-0.237 

0.144 

0.147 

-0.078 

Tfr£t 

54 aa 

0.307 0.000 

0.000 

0.605 

-0.060 -0.311 

-0.012 

1.196 

0.024 

0.235 

0.151 

-0.338 

0.079 

-1.580 

-0.145 

0.165 



• 

0.000 

0.000 

- 

- 


- 

- 

- 

- 

- 

- 

0.000 

-0.085 

Tfb>ft 

201 aa 

0.041 0.000 

0.000 

0.011 

0.000 0.000 

0.000 

-0.041 

0.000 

0.056 

0.003 

-0.055 

0.000 

0.017 

0.000 

-0.017 


273 aa 

0,084 -0.275 

0.000 

-1.724 

-0.175 -0.271 

-0.193 

0.150 

-0.308 

0.134 

0.045 

0.048 

-0.284 

0.698 

- 

- 


65 aa 

0.279 0.010 

0.000 

- 

0.010 -0.284 

0.022 

1.120 

0.022 

0.288 

0.180 

0.026 

-O.Oll 

-1.267 

-0.236 

0.151 
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TIMS teM> 


TIW (M«> 

Figure 6«: Plots of input, output, and internal variables for the 

"normal" run of Task #3. Variables plotted are the same 
as for Figure 4. The input signal is a 250 ms pulse of 
magnitude 60 N-m for the flexor. An external load of 12 N-m 
(not graphed) exists throughout the movement. The steady- 
state input signals take this fact, plus the static torque- 
angle relations, into account. Notice the large torque lost 
to muscle viscosity. 



Fn-«, Fi«-b 



TINE (MC> 

Figure 6b ; Sensitivity Trajectories 
Sensitivity Graphs (righ 
"EXT load", an external 
value of 12 N-m. Variab 
same as in Figure 5. 


















TIME (MC> 


Figure 6d, , Sensitivity Trajectories (left panel) and Sensitivity Graphs 

(right panel) for an initial velocity of -2 rad/sec (nominal). 
Axes ranges and behaviors plotted are the same as in Figure 5. 
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Task #4: Maximal Isometric Contraction: This Is a fundamentally 
different type of contraction than those considered above. A maximal 
contraction is resisted by an equal and opposite external force. An 
interesting question is whether or not a different set of parameters become 
most important. The sensitivity matrix of Table 5 shows some sharp 
differences both in what parameters and in what behaviors are most 
important. First of all, since there is no movement, there is no change in 
overall position, velocity or acceleration, leaving only force propogation. 
As expected, the passive plant parameters Jp, Bp and Kp have no effect. We 
now see that the series elasticity, as well as Hill'’s constants and the time 
constant of activation/ deactivation, play a role. Also, because, for this 
task of contraction, the hypothetical internal muscle node is shortening, 
the lengthening muscle parameters like Fm-fv also do not play a role until 
after 500 msec. 

Inspection of the sensitivity matrix indicates that only the viscous 
torque behavior seems to be effected. This result is misleading, as can be 
seen from the sensitivity trajectory for Ks2-f (Figure 7). This example 
shows one of the major limitations of sensitivity coefficients and graphs: 
in the process of extracting useful information, other information is 
invariably lost, and, furthermore, the information extracted can sometimes 
be misleading. It turns out that these observations are pairticularlly true 
for the series elastic parameters, which tend to primarily effect movement 
with high frequency components, such as movement initiation or voluntary or 
involuntary limb oscillation. Thus, subtle information on trajectory shapes 
can be lost when sensitivity coefficients and graphs are restricted to 
behaviors such as peak output values and the corresponding time of peak 
values. Although such phenomena are difficult to define by the "behaviors" 
presented here, they are possibly discribable by other types of behavior 
definitions, such as oscillation frequency. In any case, sensitivity 
trajectories must be plotted. 

TABLE 5: 


SENSITIVITY MATBIX POE TASK: "Maxiical Voluntary laoitetrlc Contraction": 




Jp 

Bp 

Kp 

Ks1-f 

Kb2-T 

Af-f 

Bh-f 

Pmfv 

Tal 

Ta2 

Pm-f: 

60.0 N-m 

0 

0 

0 

0.000 

0.000 

0.000 

0.000 

0 

- 0.002 

0.000 

Pt-f: 

15.2 N-m 

0 

0 

0 

- 0.184 

- 0.727 

0.278 

- 0.727 

0 

- 0.610 

-0.135 

Tfm-f: 

500 ms 

0 

0 

0 

0.000 

-0.001 

0.000 

-0.001 

0 

0.289 

0.000 

Tfb-f: 

43 tts 

0 

0 

0 

-0.017 

- 0.360 

0.132 

- 0.360 

0 

0.334 

0.292 


Task #5: A Simple External Load, With No Neural Control: Here we 
have a steady neural input signal of 6 N-m for both the flexor and extensor 
groups. An external load of 12 N-m is applied for 200 msec, and no effort 
is made to resist this disturbance via neural feedback. Inspection of the 
resulting sensitivity matrix (Table 6) shows again that the relative 
sensitivity of the various parameters and behaviors is a function of the 
task. 
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Figure 7; Sensitivity Trajectory for parameter Ks2, a series elastic fit 
parameter (left panel). Small dashed is for 1/5 and 5 times 
nominal - note lack of symmetry. Sensitivity Graphs (right 
panel) are for the same parameter. Ranges are from 1/10 to 10, 
in logarithmic units, for both the parameter (abscissa) and the 
behaviors (ordinate). Note the lack of information for Fm-f. 


TABLE 6 ; 

SENSITIVITY ANALYSIS FOR: External Load with no Active Resistance 



NOMINAL 

Jp 

Bp 

Kp 

Kp2 

Ksl-f 

Ks2-f 

Af-f 

Bh-f 

nafv 

Tal 

Magn: 

49.4 deg 

-0.172 

-0.154 

-0.142 

-0.093 

-0.015 

-0.032 

-0.164 

0.394 

-0.375 

-0.010 

Vmax: 

267.5 d/s 

-0.181 

-0.176 

-0.147 

-0.024 

-0.008 

-0.018 

-0.188 

0.437 

-0.405 

-0.013 

Anax: 

5730 d/s/s 

-1.000 

0.000 

0.000 

0.000 

0.000 

0.000 

0.001 

0.000 

0.001 

0.000 

Amin: 

-5742 d/s/s 

-1.108 

0.018 

0.078 

0.816 

-0.004 

-0.010 

0.018 

-0.032 

0.019 

0.012 

Tmagn: 

266 ms 

0.277 

-0.030 

-0.101 

0.000 

0.014 

0.022 

-0.049 

0.100 

-0.095 

-0.008 

Tvmax: 

196 ms 

0.468 

-0.131 

-0.275 

-0.223 

0.033 

0.076 

-0,158 

0.335 

-0.264 

-0.057 

Tamax: 

1 ms 

0.000 

0.000 

0.000 

0.000 

0.000 

0.000 

0.000 

0.000 

0.001 

0.000 

Tamin: 

200 ms 

0.004 

0.000 

0.000 

0.003 

0.001 

0.000 

0.000 

0.000 

0.000 

0.000 

Fm-f: 

6.00 N-m 

0.001 

0.000 

0.030 

0.115 

0.000 

0.000 

0.005 

-0.033 

0.041 

0.000 

Fm-e: 

7.47 N-m 

-0.018 

-0.006 

-0.030 

-0.000 

-0.001 

-0.003 

-0.006 

0.017 

0.183 

0.003 

Fb-f: 

2.87 N-m 

-0.135 

-0.094 

-0.062 

-0.004 

-0.003 

-0.008 

0.165 

-0.402 

-0.234 

0.011 

Fb-e: 

1.51 N-m 

-0.064 

-0.039 

-0.023 

0.385 

-0.004 

-0.009 

-0.041 

0.101 

0.711 

0.006 

Fm-fe; 

4.38 N-m 

-0.110 

-0.074 

-0.048 

-0.002 

-0.003 

-0.008 

0.096 

-0.227 

0.173 

0.009 

Tfm-f : 

2 ms 

- 

0.000 

— 


0.000 

0.000 

_ 

- 

- 

0.000 

Tfm-e: 

101 ms 

0.516 

-0.014 

-0.051 

-0.007 

0.029 

0.065 

-0.051 

0.115 

0.098 

0.107 

Tfb-f: 

154 ms 

0.518 

-0.071 

-0.182 

-0.048 

0.019 

0.038 

-0.087 

0.226 

-0.136 

0.047 

Tfb-e; 

135 ms 

0.591 

-0.043 

-0.133 

0.418 

0.032 

0.077 

-0.072 

0.165 

-0.858 

0.069 

T£m-fe: 

150 ms 

0.530 

-0.064 

-0.169 

-0.045 

0.019 

0.044 

-0.074 

0.196 

-0.139 

0.048 


In Figure 8 the sensitivity trajectories and graphs are obtained for an 
interesting parameter that is not well understood. This is "Fm-fv", a 

torque— velocity parameter, discussed earlier, that influences only 
lengthening muscle (by giving the lenthenlng muscle torque eccentric torque 
saturation value, nominally 30% above isometric for any given activation 
level). The fact that this parameter is significant shows that, for this 
task, the constitutive relation used for the lengthening extensor muscle 
group is important. 
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Figure 8: Sensitivity Trajectories (left panel) and Sensitivity 

Graphs (right panel) for parameter Fm-fv, a torque- 
velocity parameter for lengthening muscle. This parameter 
gives the percentage above isometric that a muscle can 
contract when lengthening at medium-to-high velocity, 
i.e. the saturating torque for lengthening muscle. 
Trajectory ranges, axes labeling, and axes ranges are 
as in Figure 5. 




DISCUSSION; 


The five simple tasks considered above barely scratch the surface of 
all the potential movements of the elbow joint. These movements seemed 
representative of the range of possibilities. Furthermore, the neural 
Inputs were purposely kept simple for the sake of clarity of presentation. 

In reality, however, neural inputs - as well as system output, are more 
complex. Over the last few decades a large amount of research has been 
done on upper limb movement, including elbow movements performed both in 
Isolation and in conjunction with movements of other joints. Much of this 
data is for atheletic performance. Components of such data can be simulated 
using JAMM. This is possible because the model contains all of the 
parameters that are needed to describe all basic nonlinear muscle 
properitles - the model is purposely constructed to be able to simulate a 
the full range of types of tasks seen in the literature. 

By using sensitivity analysis the parameters of primary importance can 
be determined for any particular task. This provides insight into the 
movement task under analysis. It also suggests ways for task-specific model 
reduction, if desired. 

Another area of Interest is the sensitivity analysis protocol. Once a 
task for analysis is chosen, the following steps were found to be represent 
an optimal protocol: First, a sensitivity aatrix is constructed, the size 
of which depends on the parameters and behaviors of Interest. The 
coefficients in the matrix are usually best found using the "logarithmic” 
method. This gives one a global view of model performance and furthermore 
guides one to the areas of interest for more detailed work. Second, 
sensitivity trajectories are used to help visualize the effect of a given 
parameter. Sensitivity Trajectories are also a good way to scan for 
problems in coefficient values. Finally, sensitivity graphs are of 
considerable use in getting a feel for the model behavior as the parameter 
is varied over a wider operating range. Such Information often descibed the 
potential tolerable range of the parameter and also the linearity of the 
change in behavior with change in parameter. As such it also can show the 
extent of the sensitivity matrlx°s sensitivity to range used for coefficient 
determination. 

There are a couple of observations worth noting: First, all three 
methods have weaknesses. Second, one of the main advantages of using such a 
variety of methods is that the weaknesses of any particular method are 
exposed by the other methods. Consequently, each method gains strength when 
combined with the others. Third, the work presented here was only for elbow 
flexion-extension. These methods are currently being found equally valuable 
for knee flexion-extension, wrist flexion-extension, and eye, head and wrist 
rotation. 
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CONCLUSION! 


A number of sensitivity analysis tools have been used on a highly 
developed model of elbow flexion-extension. It has been found that maximum 
insight into both model performance and parameter sensitivity appears to 
require a systematic protocol that employs a variety of sensitivity tools. 
Each of the methods is strengthened when used in conjunction with the other 
sensitivity tools. Furthermore, it has been found that the relative 
sensitivity of the model parameters is a function of the task being studied. 

Finally, it is suggested that sensitivity analysis should be the cornerstone 
for task-specific model reduction. 
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